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When  two  different  substances  are  joined,  material  (lows  across  the  interface  (sometimes  al¬ 
most  imperceptibly)  until  the  chemical  potentials  of  the  component  species  arc  equalized. 
When  the  substances  are  solid  or  liquid  and  some  of  the  chemical  species  are  charged,  then 
the  interface  develops  a  net  electrical  polarization  due  to  the  formation  of  an  electric  double 
layer.  The  main  goal  of  this  program  of  study  is  to  give  a  molecular  basis  for  understanding 
the  structure  and  dynamics  of  electric  double  layers  at  charged  metal-aqueous  electrolyte 
interface.  The  aim  is  to  unify  current  separate  descriptions  of  surface  adsorption  and  solution 
behavior  and,  ultimately,  to  include  a  detailed  treatment  or  the  surface  crystalography  and 
electronic  properties  of  the  metal. 
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INTRODUCTION 

When  two  different  substances  are  joined  material  flows  across  the  interface  (sometimes 
almost  imperceptibly)  until  the  chemical  potentials  of  the  component  species  are  equal¬ 
ized.  When  the  sustances  are  solid  or  liquid  and  some  of  the  chemical  species  are  charged 
then  the  interface  developes  a  net  electrical  polarization  due  to  the  formation  of  an  electric 
double  layer.  The  existence  of  electric  double  layers  was  first  recognized  by  von 
Helmholtz  1  who  studied  them  in  the  last  century.  In  many  chemical  and  biological  sys¬ 
tems  the  electric  double  layer  exerts  a  profound  effect  on  function.  For  example  in 
aqueous  electrolyte  solution  the  electric  field  of  a  charged  object  (electrode  surface  or  an 
ion)  is  completely  shielded  by  the  movement  of  ions  of  opposite  sign  toward  the  surface 
until  charge  balance  is  achieved.  The  distribution  of  ions  around  charged  objects  is  de¬ 
scribed  simply  by  the  classical  theories  of  Gouy2"4  and  Chapman5  for  flat  surfaces  and 
by  the  theory  of  Debye  and  Hiickel6  for  spherical  ions. 

The  main  goal  of  this  program  of  study  is  to  give  a  molecular  basis  for  understanding  the 
structure  and  dynamics  of  electric  double  layers  at  charged  metal-aqueous  electrolyte 
interface.  The  aim  is  to  unify  current  separate  descriptions  of  surface  adsorption  and 
solution  behavior,  and  ultimately  to  include  a  detailed  treatment  of  the  surface 
crystalography  and  electronic  properties  of  the  metal.  A  key  element  in  this  effort  is  the 
correct  treatment  of  electrostatic  interactions  among  ions,  polar  neutrals  and  their  images 
in  the  metal.  In  our  work  this  is  achieved  with  the  use  of  the  fast  multipole  method  of 
Grcengard  and  Rokhlin7'10.  The  fast  multipole  method  was  specifically  designed  for  ef¬ 
ficient  computation  of  long  range  coulomb  interactions.  In  our  simulations  we  have  used 
it  to  in  die  calculation  of  all  electrical  forces,  including  direct  space  and  image  inter¬ 
actions. 


Because  we  calculate  the  electrostatic  interactions  accurately,  without  the  use  of  finite 
cut-offs,  this  work  also  has  important  implications  concerning  the  way  in  which  ions  and 
water  distributions  in  biological11  and  clay  suspensions12  systems  should  be  calculated. 

Four  previous  publications  summarize  the  work  completed  prior  to  that  described  in  this 
report13"16.  Thmr  '  systematic  computer  simulations  we  have  shown  that  a  relatively 
aiiupie  model  suffices  to  describe  the  adsorption  of  halide  anions  and  alkali  metal  cations 
on  neutral  and  charged  metal  surfaces.  These  calculations  show  qualitatively  many  of  the 
features  known  to  occur  experimentally  at  electrochemical  interfaces  in  the 
thermodynamically  stable  region  of  the  electric  double  layer.  Most  notable  are  the  exist¬ 
ence  of:  highly  oriented  water  layer  next  to  the  charged  metal,  contact  adsorption  layer 
of  the  larger  ions,  and  a  diffuse  region  of  strongly  hydrated  ions.  The  contact  adsorbed 
ions  comprise  the  compact  layer  in  electrochemical  interfaces.  In  our  model,  contact 
adsorbed  ions  are  physisorbed  because  there  is  no  provision  to  describe  covalent  bonds. 
Typical  electric  fields  found  in  our  calculations  correspond  to  5x1 09  V/m,  and  these  arise 
from  image  charge  densities  on  the  electrode  of  about  0.1e/nm2.  One  by  product  of  our 
simulations  is  to  point  up  some  of  the  limitations  in  current  electrochemical  concepts. 
Consider  for  example  the  outer  Helmholtz  plane  which  defines  the  position  of  closest 
approach  of  hydrated  ions  to  the  elctrode  surface.  This  plane  is  shown  in  numerous 
textbook  and  review  article  illustrations17*  18.  Our  simulations  show  that  the  plane  really 
corresponds  to  a  narrow  zone  where  ions  of  the  diffuse  layer  ions  penetrate  with  increasing 
difficulty  the  closer  they  are  to  the  electrode  surface.  In  this  zone  some  mixing  with 
contact  adsorbed  species  can  occur  which  is  not  possible  in  the  old  picture. 

In  the  last  fifteen  years  there  have  been  many  simulations  of  water  and  electrolyte  sol¬ 
utions  near  surfaces.  Some  of  these  studies  have  contributed  greatly  to  our  understanding 
of  electrochemical  interfaces.  For  completeness  some  of  this  work  is  summarized  here. 
Films  of  pure  water  between  uncharged  dielectric  walls13*  19"21,  and  charged  dielectric 
walls13*  14*  22*  23.  Some  of  this  work  is  noteworthy  because  of  a  predicted  phase  tran¬ 
sition22*  23.  There  have  been  numerous  calculations  reported  for  uncharged  metal24-27 
walls28"32,  including  one  for  jellium30  and  several  for  corrugated  platinum 
surfaces27-29*  31*  32  predicting  that  water  adsorbs  weakly  at  top  sites  with  oxygen  down 
on  Pt(lll)  and  Pt(100).  There  have  also  been  some  calculations  for  electrolyte  solutions 
between  uncharged  and  charged  dielectric  walls13*  14*  33*  34  emphasizing  spatial  distrib¬ 
utions  and  hydration  shell  structure.  There  have  been  studies  for  electrolytes  between 
uncharged  metal  walls27*  35*  36.  The  work  of  Rose  and  Benjamin36  is  particularly  inter¬ 
esting  because  umbrella  sampling  was  used  to  calculate  the  free  energy  of  adsorption. 
Finally  we  mention  the  studies  of  water  between  charged  metal  walls37,  and  electrolytes 
between  charged  metal  walls37.  In  much  but  not  all  of  the  wotk  just  summarized,  the  long 
range  coulomb  interactions  were  treated  in  an  approximate  way.  The  commonest  ap¬ 
proximation  was  to  cut  off  all  interactions  beyond  a  certain  radius  like  0.10  nm.  Some 
workers  used  the  Ewald  method  or  a  planar  modification  of  the  method  to  compute  the 
sum  of  long  range  fields  correctly.  Curiously  many  of  the  other  summation  methods,  like 
the  planewise  method  of  de  Wette  38  seem  not  to  have  been  used  at  all.  The  reason  why 
it  is  important  to  calculate  the  long  range  fields  accurately  is  to  capture  the  macroscopic 
part.  The  dipole  component  of  the  field  is  conditionally  convergent  and  the  sum  must  be 
performed  in  a  manner  consistent  with  the  physical  boundaries.  This  is  a  rather  old 
problem  which  has  a  partial  mapping  onto  the  problem  of  calculating  the  macroscopic 
electromagnetic  field  inside  a  sample  of  arbitrary  shape.  Space  and  time  (!)  prevent  us 
from  pursuing  this  connection  here. 


CAPILLARY  THERMODYNAMICS 


Our  present  understanding  of  electric  double  layers  and  their  importance  in  surface 
electrochemistry  has  its  origin  in  the  dropping  mercury  electrode  experiments  of 
Grahame39’  40.  The  basic  experiment  consists  of  measuring  the  radius  of  a  sessile  mer¬ 
cury  drop  as  a  function  of  electrode  potential  and  electrolyte  composition.  One  can  then 
use  the  electrocapillary  equations  of  Lippmann41  to  relate  radius  of  the  drop  to  the 
interfacial  surface  tension  and  so  in  turn  to  the  charge  on  the  electrode.  Thermodynamic 
arguments  are  then  used  deduce  the  surface  excess  concentration  of  adsorbed  cations, 
anions  and  neutral  organics  on  the  electrode  surface17*  18*  39*  40, 42*  43. 

Thermodynamics  of  Adsorption 


electrocapillary 
thermodynamics 
Lipmann  eqn. 


non-specific 
adsorption 
PFj,  F*.  CIO4,  ... 
Li+,  Na+, ... 
neutrals 


contact  adsorption 
I*.  Br,  Cl',  SCN',  N3, ... 
NR4.  Cs, ... 


Figure  L  Schematic  diagram  depicting  thermodynamics  of  adsorption:  (a)  plot  of  surface  tension  y  vs 
electrode  potential  V  for  different  electrolyte  compositions,  (b)  surface  concentration  T  s  vs  potential  V  for 
neutrals  and  strongly  hydrated  ions,  (c)  Ts  vs  potential  V  for  contact  adsorbing  ions. 


Figure  1  summarizes  the  important  aspects  of  thermodynamic  studies.  For  example  at  the 
top  of  Figure  1  the  surface  tension  y  is  depicted  versus  electrode  potential  for  aqueous 
NaF  at  different  concentrations.  The  Lippmann  equation 


[1] 


can  be  used  to  determine  the  charge  on  the  mercury  drop  directly,  and  in  turn  the  differ¬ 
ential  capacitance  can  be  calculated.  Contact  adsorbed  ions  give  a  large  contribution  to 
the  capacitance  compared  to  ions  in  the  diffuse  layer,  and  this  has  be  used  to  separate  their 


contributions  in  mixed  solutions.  The  middle  Figure  1  shows  highly  schematic  curves  for 
the  adsorption  of  neutral  organics  and  strongly  hydrated  ions  on  noble  metals.  The  or¬ 
ganics  adsorb  most  strongly  when  the  charge  on  the  electrode  is  small  because  the  water 
layer  nearest  the  metal  is  then  least  strongly  bound  and  more  easily  displaced  by  the  or¬ 
ganic  adsorbate.  Hie  final  part  Figure  lc  shows  the  adsorption  isotherms  for  larger  ions, 
in  particular  it  depicts  behaviour  like  that  observed  for  iodide  and  the  pseudohalides  CN~ 
and  thiocyanate  SCN\  For  more  discussion  see  the  review  paper  by  Anson44.  For  solid 
electrodes  less  direct  methods  are  used  to  obtain  the  relevant  data  for  a  thermodynamic 
analysis.  If  the  differential  capacitance45,  46 


C(V)  = 


[2] 


can  be  measured  in  a  way  that  inner  layer  and  diffuse  layer  contributions  can  be  separated 
then  a  variety  of  integration  techniques  can  be  used  to  determine  the  charge  on  the 
electrode.  From  these  quantities  the  other  important  thermodynamic  variables  such  as  the 
surface  concentrations  Ts  (depicted  schematically  in  middle  and  bottom  of  Figure  1)  can 
be  obtained  using  the  thermodynamic  analysis  starting  from  Eqn  (1). 


THE  TRADITIONAL  MODEL 

On  the  basis  of  experimental  methods  like  differential  capacitance,  chronocoulombmetry, 
ellipsometry,  and  UV-visible  spectroscopy,  wielded  with  consummate  skill  a  detailed 
picture  of  the  electric  double  layer  adjacent  to  metal  surfaces  has  been  devised18-  47 '. 
We  will  refer  to  this  picture  as  the  'traditional  model'.  The  main  features  of  this  model 
are  shown  schematically  in  Figure  2  for  that  part  of  the  double  layer  close  to  the  electrode 
surface.  Drawings  similar  to  this  one  can  be  found  in  many  textbooks  and  review  articles 
on  interfacial  electrochemistry48*50.  The  metal  is  flat  and  carries  charge  (negative  shown), 
the  aqueous  subphase  is  divided  into  two  parts,  called  diffuse  and  compact  regions.  Ad¬ 
ditionally  next  to  the  charged  surface  there  is  a  highly  oriented  layer  of  water  shown  in 
Figure  2  with  protons  oriented  towards  the  surface.  The  anion  is  shown  adsorbed  in 
contact  (physisorbed)  and  the  strongly  hydrated  cations  are  shown  no  closer  than  two 
water  molecules.  The  inner  Helmholtz  plane  (EHP)  is  defined  as  the  plane  through  the 
nuclei  of  the  contact  adsorbed  anions  and  a  similar  plane  through  cations  at  their  distance 
of  closest  approach  is  called  the  outer  Helmholtz  plane  (OHP).  Beyond  the  OHP  the 
distribution  of  ions  is  assumed  to  be  described  by  the  Gouy-Chapman  theory,  which  in  it's 
simplest  form  assumes  the  ions  are  charged  point-like  objects  and  the  solvent  is  a 
dielectric  continuum  with  appropriate  bulk  properties.  Close  to  the  electrode  the  'tradi¬ 
tional  model'  calculates  system  properties  based  on  static  distributions  or  uses  lattice  sta¬ 
tistics.  The  diffuse  region  in  this  picture18  starts  two  solvent  molecules  from  a  flat 
electrode  surface  and  stretches  out  several  nanometers  into  the  bulk  electrolyte.  The 
electrostatics  and  ionic  distributions  in  this  diffuse  part  were  first  described  by  the  Gouy- 
Chapman  theory2*5  which  predates  even  the  Debye-Hiickel6  model  of  ionic  atmospheres 
in  bulk  electrolytes.  The  ions  in  the  diffuse  layer  screen  the  net  charge  of  metal  and  any 
ions  occupying  in  the  compact  part  of  the  double  layer.  Traditionally  the  structure  of  the 
compact  region  is  thought  of  as  being  rather  static  and  resembling  a  parallel  plate  capacitor 
with  a  gap  of  atomic  separations  (0.1  -  0.2  nm).  The  flat  surface  model  dates  from  times 
before  the  ability  to  make  useable  single  crystal  electrode  surfaces,  the  advent  of 
synchrotron  sources  and  facilities  to  do  surface  X-ray  analysis.  Even  so  quite  a  consid- 


erable  effort  has  been  directed  toward  developing  better  model  of  the  metal  side.  In  this 
regard  the  work  of  Halley  et  al51*53  and  Schmickler54  is  particularly  significant.  These 
groups  have  developed  theories  based  on  jellium  electrode  models  of  the  charged  metal 
surface  that  attempts  to  capture  important  physics  causing  features  in  measured 
capacitance  vs.  potential  curves.  Focussing  more  on  the  electrolyte  subphase  Henderson 
and  coworkers55  have  developed  a  correlation  function  approach  and  an  analysis  that 
shows  that  in  their  very  reasonable  model  there  is  no  sharp  division  between  the  diffuse 
and  inner  Helmholtz  planes.  This  is  in  contrast  to  the  ideas  conveyed  by  the  pictures 
found  in  many  textbooks  and  Figure  2. 


Figure  2.  Schematic  diagram  depicting  the  'traditional  model'  of  the  electric  double  layer  found  in  many 
review  articles  and  introductory  texts  on  electrochemistry. 

It  must  be  pointed  out  that  many  of  molecular  dynamics  studies  performed  to  date  have 
not  recognized  the  importance  of  solving  the  electrostatics  problem  accurately.  Intense 
long  range  fields  exist  at  the  interface  because  of  the  large  dipole  normal  to  the  surface 
formed  by  the  ions  and  their  electrical  images  are  all  in  phase.  This  problem  spills  over 
into  other  areas  where  electrostatics  is  important,  for  example  most  studies  of  water  around 
biological  objects  are  wrong  in  the  way  in  which  the  dielectric  polarization  of  water  is 
calculated  because  long  range  correlations  are  lost  when  interactions  are  cut  off  at  a  finite 
distance.  Workers  in  this  field  often  resort  to  using  distance  dependent  dielectric  functions 
to  shield  ionized  groups  attached  to  the  object's  surface. 

The  experimental  picture  of  electrochemical  surfaces  is  currently  undergoing  rapid  change 
due  to  numerous  advanced  in  situ  and  ex  situ  UHV  surface  science  probes  of  the 
electrochemical  interface.  The  new  synchrotron  based  X-ray  surface  crystal  structure 
probes  like  grazing  incidence  X-ray  scattering  (GIXS),  X-ray  standing  wave  (XSW),  sur¬ 
face  extended  X-ray  absorption  fine  structure  (SEXAFS)  techniques  have  allowed  surface 
geometries  to  be  measured  for  the  first  time56-5^.  Recent  studies  by  Hubbard  and 


coworkers46  reveal  changes  in  charge  state  of  the  ion  as  the  surface  concentration  is  in¬ 
creased,  changes  that  have  analogies  in  UHV  surface  science  studies  of  die  metalization 
of  semiconductors  These  studies  show  contact  adsorption  to  be  a  complex  process  that  can 
evolve  into  chemisorption  at  high  coverages  even  if  the  initial  step  is  physisorption.  Local 
probes  like  scanning  tunneling  microscopy  (STM)  and  atomic  force  microscopy  (AFM) 
give  images  with  local  atomic  scale  features.  The  next  few  years  will  likely  see  major 
revisions  in  our  experimental  understanding  of  the  less  dynamic  part  of  the 
electrochemical  double  layers. 


MODEL  FOR  THE  IMMERSED  ELECTRODE 


Consider  a  system  consisting  of  two  electrodes  immersed  in  aqueous  electrolyte  solution 
as  shown  schematically  in  Figure  3.  Reading  from  left  to  right  there  are  three  regions: 
the  anode  on  the  left,  the  bulk  region  at  the  center,  and  the  cathode  region  on  the  right. 
If  the  electrodes  are  uncharged  then  in  all  three  regions  the  electrolyte  phase  would  be 
electrically  neutral.  Now  when  the  externally  applied  (battery)  potential  is  altered  so  that 
the  electrodes  become  charged  (the  case  shown  schematically  in  Figure  3)  the  electrolyte 
responses  by  screening  the  electric  field  and  the  three  regions  acquire  different  net  charge. 
One  with  excess  anions  (left)  screens  the  positive  charge  on  the  left  electrode,  a  bulk  re¬ 
gion  in  which  the  electric  field  is  zero,  and  one  with  an  excess  of  cations  (right).  When 
the  charge  on  each  electrode  is  included  with  the  adjacent  'zone'  of  electrolyte  the  net 
charge  in  each  region  is  zero. 

These  simple  considerations  suggest  we  can  tiy  to  model  an  immersed  electrode  with  it's 
adjacent  screening  region  and  do  not  have  to  model  the  whole  cell  as  was  done  in  two 
previous  publications13*  14.  This  approach  is  useful  because  it  reduces  the  number  of 
water  molecules  in  the  calculation,  however  it  imposes  a  constraint  in  the  form  of  charge 
neutrality  and  requires  us  to  choose  a  'good'  boundary  to  separate  the  bulk  region  and 
immersed  region.  We  can  test  the  boundary  by  scaling  the  size  of  the  system  to  capture 
some  bulk-like  behavior.  Scaling  by  a  factor  of  ten  suggests  that  four  layers  of  water  is 
a  minimum  sized  system. 

Our  immersed  electrode  model  therefore  consists  of  a  layer  of  electrolyte  between  two 
walls.  The  wall  on  the  left  carries  no  charge  it  is  simply  a  restraining  wall  and  ideally 
should  allow  a  continuous  transition  to  the  bulk  electrolyte  region.  The  complete  system 
of  electroltye  and  electrode  (always  chosen  to  be  on  the  right  hand  side  in  our  calculations) 
is  neutral,  unless  we  deliberately  choose  to  stress  the  system  with  an  uncompensated 
charge  on  the  right  hand  metal  electrode.  In  practice  we  have  performed  both  types  of 
calculation  for  reasons  described  in  more  detail  later.  Later  we  will  show  an  important 
result  that  water  behaves  in  an  uncompensated  field  as  if  the  excess  ionic  charge  were 
higher. 

Since  there  is  integer  charge  in  the  aqueous  phase  qAq,  and  the  image  charge  on  the  metal 
qim  satisfies  the  equation  qim  +  <?*,  =  0.  In  our  calculations  we  also  consider  cases  (to  be 
discussed  at  length  later)  where  there  is  additional  charge  quc  on  the  electrode  with  the 
restriction  that  qr  =  qu*  +  q*c  =  ne  where  n  =  0, ±  1,±2, ... . 

In  Figure  3  if  one  thinks  of  the  vertical  dash  line  as  symbolizing  the  restraining  wall  then 
in  reality  the  vertical  line  will  be  very  close  to  the  corresponding  electrode  for  a 
macroscopic  sized  cell.  Integral  electrode  charge  is  an  essential  constraint  in  this  im¬ 
mersed  electrode  model.  We  stress  again  that  the  main  advantage  of  the  model  is  that 


only  about  half  the  number  of  water  molecules  are  needed  to  simulate  a  system  with  two 
metal  electrodes.  For  the  restraining  barrier  we  choose  a  9-3  potential.  The  origin  for  this 
restraining  potential  is  at  1.862  nm  from  the  image  plane  of  the  metal.  We  refer  to  this 
restraining  potential  surface  as  the  dielectric  surface,  and  as  already  mentioned  it's  only 
function  is  to  limit  the  extent  of  the  fluid  phase  and  thereby  make  the  calculations  more 
tractable.  In  all  the  calculations  reported  here  the  simulation  cell  was  a  cube  with  edge 
1.862nm.  The  cube  was  periodically  replicated  in  the  xy  directions  parallel  to  the 
electrode  surface  plane.  Again  we  mention  that  we  have  performed  some  scaled  up  cal¬ 
culations  on  cells  with  larger  edge  length  by  up  to  a  factor  2  and  found  very  similar  fea¬ 
tures  to  those  described  for  the  smaller  cells. 


Figure  3.  Immersed  Electrode  Model.  Schematic  diagram  of  the  electrochemical  cell  showing  the  im¬ 
mersed  electrode  configuration  on  the  right  side.  Vertical  broken  lines  symbolize  the  transition  region  from 
diffuse  layer  to  bulk  electrolyte  where  the  solution  is  electrically  neutral. 


MODEL  FOR  WATER,  IONS  AND  THE  METAL  SURFACE 

In  all  the  calculations  reported  here  we  have  used  the  parameters  of  the  Stillinger60,  61 
ST2  water  model  and  the  extensive  interaction  parameter  set  for  alkali  metal  ions  and 
halide  ions  developed  by  Heinzinger  and  coworkers62.  Figure  4  shows  a  schematic  of  the 
ST2  water  model,  a  simple  ion  and  the  smoothly  truncated  Lennard-Jones  potential.  The 
ST2  water  molecule  model  consists  of  a  central  oxygen  atom  (0_ST2  or  O  for  short) 
surrounded  by  two  hydrogen  atoms  (H_ST2  or  H  for  short)  and  two  massless  point 
charges  (PC  ST2  or  PC  for  short)  in  a  rigid  tetrahedral  arrangement  (bond  angle  = 
cos'^l/yT ) ).  The  O-H  and  O-PC  bond  lengths  were  0.10  nm  and  0.08  nm  respectively. 
This  small  difference  in  bond  lengths  means  that  the  water_ST2  model  and  its  electrostatic 
image  (i.e.,  q  — »  -q  )  behave  similarly.  The  only  Lennard-Jones  'atom'  in  ST2  model  is 
the  oxygen  atom.  The  hydrogen  H  and  point  charges  PC  interact  with  their  surroundings 
(i.e.,  other  atoms  and  surfaces)  by  Coulomb  interactions  only.  Their  charges  are  qH= 
0.23570lei  and  q^  -qH.  The  O  atom  carries  no  charge.  The  alkali  metal  ion  and  halide 
ion  were  treated  as  non-polarizable  Lennard-Jones  atoms  with  central  point  mass  and 


charge.  Figure  4  shows  an  ion  schematically  next  to  the  water  model.  The  atom-atom 
interaction  parameters  are  taken  from  Heinzinger's  review28*  62. 

Next  we  describe  the  interaction  between  water  and  ions  and  the  metal  and  restraining 
wall.  The  metal  was  represented  by  two  linearly  superimposed  potentials.  Pauli  repulsion 
and  dispersive  attractive  interactions  were  modelled  by  a  9-3  potential,  and  the  interaction 
of  charges  with  the  conduction  electrons  by  a  classical  image  potential.  In  the  calculations 
described  here  the  image  plane  and  origin  plane  of  the  9-3  potential  were  coincident.  This 
was  tantamount  to  choosing  die  image  plane  and  the  nuclear  plane  of  the  metal  surface 
to  be  coincident  This  was  acceptable  in  our  scheme  because  the  Lennard-Jones  core 
parameters  a  are  all  large  and  the  'thickness'  of  the  repulsive  wall  is  also  large  (ca.  0.247 
nm). 


(a)  Water/ST2  (b)  lons/ST2  (c)  Lennard-Jones 
Model  Model  Potential 


Figure  4.  The  Ion  and  Water  Models.  Schematic  diagram  summarizing  the  key  features  of  the  model  for 
ST2  water  and  the  ions.  Water  and  ions  are  treated  as  non  polarizable  Lennard-Jones  atoms  with  embedded 
charges.  Shown  on  the  left  is  a  schematic  drawing  of  the  Lennard-Jones  potential  with  smooth  cut-off. 


CALCULATION  OF  THE  ELECTROSTATIC  INTERACTIONS 

Electrostatic  sums  are  conditionally  convergent  and  great  care  must  be  exercised  in  using 
cut-offs  and  different  boundary  conditions.  The  Coulomb  field  computation  grows  as  N2 
(N  number  of  charged  particles)  unless  special  measures  are  adopted.  We  use  the  fast 
multipole  method  (fmm)  devised  by  Greengard  and  Roklin7,  in  which  the  cpu  time  grows 
as  N.  The  cross  over  point  in  efficiency  for  direct  sum  versus  fmm  can  be  as  small  as 
N  =  1000  (about  250  ST2  water  molecules).  The  original  papers  of  Greengard  should  be 
consulted  for  details  of  this  clever  algorithm.  Some  discussion  of  the  use  of  fmm  in 
electrochemical  simulations  has  been  given  by  Glosli  and  Philpott15. 


TOTAL  POTENTIAL  ENERGY  OF  THE  SYSTEM 

The  Coulomb  interaction  between  molecules  was  represented  as  sum  of  1/r  interactions 
between  atomic  point  charges.  For  the  ST2  water  model  these  interactions  result  in  H- 
bonds  that  are  too  strong  when  the  PC  and  H  atoms  are  close.  These  interactions  were 
softened  for  small  molecular  separation  in  the  way  described  by  Stillinger60  and  Lee  et 


al19  by  introducing  a  switching  function  S  that  modifies  specified  atomic  coulomb  inter¬ 
actions  at  small  separations.  The  short  range  part  of  the  intermolecular  interaction  was 
modeled  by  Lennaid-Jones  potential  between  atom  pairs  on  each  molecule.  All 
molecule-molecule  Lennaid-Jones  type  interactions  were  cut-off  in  a  smooth  fashion  at  a 
molecular  separation  R  =  0.68  nm  by  a  truncation  function  T.  The  atoms  of  each  molecule 
also  interacted  with  the  surfaces  at  z  =  ±z„  where  za  =  0.931  nm.  Both  surfaces  were 
treated  as  flat  featureless  plates  with  a  uniform  electric  charge  density  of  a  on  the  metal 
plate  at  +  zQ  if  there  is  an  uncompensated  charge  on  the  metal,  otherwise  a  =  0  This  gave 
rise  to  a  uniform  electric  field,  E  -  AnKa,  in  the  z-direction  where  K  the  electrostatic 
coupling  constant  had  the  value  138.936  kJ.nm/(mole.e2  )  in  the  units  of  this  paper.  The 
complete  interaction  energy  U  was  given  by  the  following  formula 


where  j  and  j  were  molecular  indices,  and  a  and  P  were  atomic  indices.  The  symbol  A 
represented  the  set  of  all  atoms  of  molecular  i .  The  symbol  Ry  was  the  distance  between 
the  center  of  mass  of  molecules  i  and  j.  The  symbol  rap  was  the  distance  between  atoms 
a  and  p.  For  small  R  we  followed  the  practice  of  modifying  the  the  coulomb  energy 
between  ST2  molecules  and  ions  by  the  switching  function  S(R,  RL,  Ru )  given  by. 


S(R,  Rl,  Rfj) 


0 


( R-Rl?(3Ru-2R-Rl ) 

(.Rv-Rj 

1 


r<rl 

rl<r<ru 

Ru<R 


[4] 


The  values  of  RL  and  Ru  were  dependent  on  the  types  of  the  molecular  species  that  were 
interacting.  As  mentioned  above  the  tails  of  the  Lennard-Jones  pair  interactions  were  cut 
off  by  the  truncation  function  T.  The  form  of  T  was  given  by, 


m  = 


T 

r<rl 

rtl<r<rI 

4<r 


[5] 


The  same  truncation  function  has  been  applied  to  all  non  Coulombic  molecular  inter¬ 
actions,  with  Rl= 0.63  nm  and  R{/= 0.68  nm.  The  integers  m  and  n  controlled  the  smooth¬ 
ness  of  the  truncation  function  at  R[  and  Rl  respectively.  In  this  calculation  n  =  m  =  2 
which  insured  that  energy  has  continous  first  spatial  derivatives. 

All  the  atom-atom  and  atom-surface  interaction  parameters  are  given  in  Table  I.  For  ex¬ 
ample  we  see  that  the  (e,o)  pairs  are  (0.3164,  0.3100),  (0.1490,  0.2370)  and  (0.4080, 
0.5400)  for  0_ST2,  Li  ion  and  I  ion  respectively.  The  units  of  the  well  depth  e  are 
kJ/mole  and  the  van  der  Waals  a  is  in  nm  .  The  usual  combining  rules  were  enforced  for 
unlike  species,  namely:  Eab  =  (t/ufoa)*  and  Gab  =  1 6(gm  +  o m).  The  st2  model  switching 


function  (see  later)  interval  ends  R?l  and  R%  both  vanish  except  for  st2/st2  water  pairs, 
where  /?£”’sr2=0.20160  nm  and  RfJ7,sn=Q.3\2&7  nm. 

The  atom-surface  parameters  describing  interaction  with  nonconduction  electrons  were 
chosen  to  be  die  same  as  those  used  by  Lee  at  al19,  A=17.447xl0'6  kJ(nm)6/mole  and 
B=76.144xl0‘3  kJ(nm)3/mole  for  O,  I  ion  and  Li  ion.  The  A  and  B  parameters  for  H_st2 
and  PC_st2  were  set  to  zero.  The  potential  corresponding  to  these  parameters  describe  a 
graphite-like  surface.  Real  metals  would  have  much  larger  e  ours  were  deliberately  chosen 
to  be  small  so  as  to  permit  coulombic  interactions  to  dominate  the  physics. 


Table  I.  The  interaction  parameters  ( q ,  e,  o,  A,  B)  and  mass  (m)  for  all  atoms  used  in  the  simulations,  where 
q0  =  e  ,  £»  =  l  KJ I  mole  ,  o„  =  1  nm  ,  A>  =  17.447  x  10 ~6KJ  •  nmVmole  ,  B0  =  76.144  x  \QriKJ  •  nm3fmole  .  and 
mo  -  l AMU.  The  e  and  0  for  unlike  atom  pairs  are  formed  using  the  combination  rules  £ab  =  and 

Gab  =  (Gaa  +  0m)/ 2. 


<?/<?« 

tit* 

O/O0 

m/mo 

A/Ao 

B/B0 

0_ST2 

0.0000 

0.316 

0.310 

16.0 

1 

1 

H_ST2 

0.2357 

0.000 

0.000 

1.000 

0 

0 

PC_ST2 

-0.2357 

0.000 

0.0 

0.0 

0 

0 

Li 

1.0000 

0.149 

0.237 

6.9 

1 

1 

F 

-1.0000 

0.050 

0.400 

19.0 

1 

1 

Cl 

-  1.0000 

0.168 

0.486 

35.5 

1 

1 

Br 

-1.0000 

0.270 

0.504 

79.9 

1 

1 

I 

-1.0000 

0.408 

0.540 

129.9 

I 

I 

In  the  equations  of  motion  bond  lengths  and  angles  were  explicitly  constrained  by  a 
quaternion  formulation  of  the  rigid  body  equations  of  motion63.  The  equations  of  motion 
were  expressed  as  a  set  of  first  order  differential  equations  and  a  fourth  order  multi-step 
numer.  J  scheme  with  a  2  fs  time  step  was  used  in  the  integration.  At  each  time  step 
a  small  scaling  correction  was  made  to  the  quaternions  and  velocities  to  correct  for  global 
drift.  Also  the  global  center  of  mass  velocities  in  the  x  and  y  directions  was  set  to  zero 
at  each  time  step  by  shifting  the  molecular  translational  velocities. 

In  the  analysis  of  configurations  the  first  100  ps  were  used  to  equilibrate  the  system  and 
subsequent  configurations  were  used  in  compute  properties  like  density  profiles.  There 
were  exceptions  where  it  was  evident  that  the  system  had  not  equilibrated.  In  practice  it 
was  found  this  occurred  frequently  in  three  ion  systems,  as  in  the  case  of  a  cation  in  the 
presence  of  two  coadsorbed  iodide  ions.  All  the  simulations  were  run  to  1000  ps  or 
longer,  so  that  generally  900  or  more  configurations  one  ps  apart  were  used  to  calculate 
averages.  Typical  density  plots  were  derived  from  bining  configurations  stored  every  ps 
and  with  bin  widths  of  0.005  nm  or  larger. 


EQUATIONS  OF  MOTION 


* 


All  molecules,  in  this  study,  were  assumed  to  be  a  rigid  collection  of  point  atoms,  so  that 
all  bond  lengths  and  bond  angles  within  a  molecule  were  fixed.  To  evolve  a  collection 
of  these  molecules  a  quaternion  formulation  of  the  rigid  body  equations  of  motion  was 
used63-65  The  center  of  mass  position  (R<)  and  velocity  (V,)  was  used  to  describe  the 
translational  degrees  of  freedoms  of  molecule  i.  The  orientational  motion  of  the  molecule 
was  described  by  the  quaternion  q(  =  (qf,  q\,  q},  qf)  and  the  rotational  velocity  (atf),  as 
measured  in  the  body  flame  of  the  molecule.  The  one  exception  to  this  was  for 
monatomic  molecules,  in  which  case  the  orientational  degrees  of  freedom  were  not 
needed. 

The  discussion  of  the  equations  of  motion  begin  by  considering  the  potential  energy  U. 
From  Eqn  3,  it  can  be  seen  that  the  potential  energy  can  be  treated  as  a  scalar  function 
of  the  variables  (Ru...,  R*;  ri,...,  r„),  where  N  is  the  number  of  molecules  and  n  is  the 
number  of  atoms  in  the  entire  system. 

U  =  f/(R],...,  Rw;  r„)  []6] 

Of  course  not  all  these  variables  are  independent.  However  for  the  purposes  of  the  fol¬ 
lowing  two  definitions  they  are  treated  as  independent  variables. 

=  -  V,V.  F(  — V^.  [7] 

The  total  force  FT  and  torque  xt,  can  be  expressed  in  terms  of  fa  and  F,  as 

F.r=F,+ 2f«  *,=  X<r«-Ri)xf-  w 

ae  Aj  ae  A, 


The  translational  motion  of  the  molecule  is  described  by  the  first  order  ordinary  differen¬ 
tial  equation, 


dR{ 

nr 


dVi 
M:  — ~7~~ 

1  dt 


[9] 


For  the  rotational  motion  it  is  convenient  to  work  in  the  body  axis  of  the  molecule,  where 
the  moment  of  inertia  t  is  diagonal  and  time  independent.  It  is  useful  to  define  the  op¬ 
erator  Qi 


\  2  3 

~qi  -q i  -qi 

0  3  2 
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[10] 


and  die  body  frame  rotational  force 


©f  x 


[n] 


Using  these  quantities  above,  the  following  first  order  ordinary  differential  equations  can 
be  written  to  describe  the  rotational  motion 


dqy 

dt 


[12] 


This  set  of  equations  conserve  the  total  energy  for  time  independent  potentials  (U)  .  A 
constant  temperature  ensemble  may  be  simulated  by  introducing  a  velocity  dependent  term 
in  the  acceleration  terms  to  constrain  the  total  kinetic  energy.  The  constant  temperature 
equations  of  motion  are  written  as 
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~dT 


dtu 

dt 


=  Ff-yv, 
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Y=(^vIF[  +  ©fFt)/(2^ 

i 


This  choice  of  y  ensures  that  the  total  kinetic  energy  of  the  systems  is  constant. 


INTRODUCTION  TO  THE  RESULTS 

The  purpose  of  this  section  is  to  briefly  introduce  and  explain  some  of  the  terms  and 
concepts  used  in  the  description  and  discussion  of  the  results.  We  define  inner  layer,  inner 
surface  field  and  the  external  or  as  we  will  also  refer  to  it  the  'uncompensated'  field. 

In  our  calculations  the  inner  layer  is  defined  to  be  all  ions  and  molecules  in  contact  with 
the  electrode.  Recall  that  the  surface  is  represented  by  a  9-3  potential  that  acts  on  the 
center  of  the  molecules  and  ions.  This  means  that  in  our  model  the  inner  most  layer,  i.e., 
the  first  layer,  has  a  simpler  structure  because  steric  effects  due  to  dimension  perpendicular 
to  the  surface  are  suppressed.  In  a  more  realistic  model  Pauli  repulsion  would  ensure  that 
the  center  of  water  was  closer  to  the  surface  that  the  center  of  an  iodide  ion.  One  conse¬ 
quence  of  our  model  is  that  some  hydrophobic  effects  at  the  surface  are  accentuated,  and 
some  steric  effects  based  on  volume  are  decreased  in  importance. 

In  equilibrium  states  the  charge  on  the  flat  electrode  is  just  the  net  image  charge,  which 
is  equal  in  magnitude  but  opposite  in  sign  to  the  sum  of  the  charges  on  the  ions.  In  this 
case  there  is  no  field  in  the  bulk  because  the  fields  from  the  ions  and  their  images  cancel. 
When  the  charge  on  die  electrode  is  only  the  image  charge,  then  the  potential  drop  occurs 


between  the  metal  surface  and  the  ions.  There  is  essentially  no  electric  field  beyond  the 
ions.  We  call  the  field  generated  by  the  ions  and  their  images  the  inner  surface  field. 
Since  die  water  molecules  have  large  permanent  dipole  moments  those  molecules  between 
the  ions  and  the  metal  will  try  to  align  with  the  field. 

It  is  useful  from  the  theoretical  point  of  view  to  stress  the  system  by  placing  extra  un¬ 
compensated  charge  on  the  electrode.  In  this  case  there  is  a  field  across  the  sample,  just 
as  occurs  when  we  have  water  without  ions  between  charged  plates.  We  refer  to  this  field 
as  the  external,  applied  or  uncompensated  electric  field.  In  general  then  we  can  think  of 
the  total  charge  on  the  electrode  as  the  sum  of  the  image  charge  and  the  extra  uncom¬ 
pensated  charge.  In  passing  we  emphasize  again  that  in  electrochemical  systems  it  is  not 
possible  to  have  uncompensated  charge  (i.e.,  unshielded)  on  an  electrode  immersed  in 
electrolyte.  However,  it  is  possible  in  principle  to  place  uncompensated  charge  on  an 
emersed  electrode,  so  that  there  are  physically  realizable  states  resembling  those  of  the 
immersed  model  with  uncompensated  charge.  However  the  uncompensated  charge  den¬ 
sities  reachable  in  practice  are  about  100  times  smaller  than  the  ones  used  in  this  paper. 
In  our  calculations  with  uncompensated  charge  we  have  deliberately  chosen  charge  den¬ 
sities  equal  in  magnitude  to  the  image  charge  density  in  order  to  better  understand  the 
effect  of  the  surface  field  on  the  inner  layer  of  water. 

Finally  we  state  again  that  the  inner  surface  field  is  localized  close  to  the  metal,  in  contrast 
to  the  field  of  uncompensated  charge  which  extends  across  the  entire  sample.  At  the 
surface  the  total  surface  field  consists  of  the  inner  field  and  any  applied  field  due  to  un¬ 
compensated  charge. 


POSITIVELY  CHARGED  METAL  IN  WATER  WITHOUT  IONS 

Water  is  the  common  primary  chemical  species  in  all  the  calculations  described  in  this 
paper  and  it  is  important  to  know  its  behaviour  in  electric  fields  without  the  additional 
structure  changing  effects  of  ions.  In  this  simulation  158  ST2  model  water  molecules 
were  confined  between  the  metal  surface  at  z  =  0.931  nm  and  the  'dielectric'  wall  at  z 
=  -zQ.  Gap  between  surfaces  is  Az  =  1.862  nm  =  L  the  edge  length  of  the  simulation  cell. 
The  water  film  thickness  corresponds  to  about  4.5  layers  of  water.  For  comparison  we 
note  that  in  the  calculations  of  Lee,  McCammon  and  Rossky19  there  were  216  water  ST2 
molecules  equivalent  to  about  six  layers.  The  image  plane  of  the  metal  was  at  z  -  0.931 
nm  the  right  hand  surface.  Left  confining  boundary  is  a  dielectric  surface  (with  dielectric 
constant  e=l)  with  no  image  field.  The  repulsive  part  of  the  9-3  potentials  on  both  sides 
of  die  box  began  at  lzl=0.682  nm. 

There  are  no  ions  in  this  sample.  The  electric  field  across  the  system  comes  from  (un¬ 
compensated)  positive  charge  density  of  <?r  =  <?«:  =  + 1  on  plate  with  area  L2,  or  equiv- 


alently  an  charge  density  of  0.29e/(nm)2.  The  electric  Held  in  vacuum  due  to  this  charge 
is  about  5.2  GV/m  .  Density  profiles  in  the  z  direction  were  obtained  by  averaging  900 
configurations  over  the  xy  plane.  These  profiles  are  shown  in  Figure  5.  The  total  surface 
field  is  just  the  applied  field.  Since  the  field  attracts  negative  charge  to  the  metal  the  water 
orients  with  point  charges  PC_ST2  directed  towards  the  metal.  In  comparing  Figure  5 
profiles  with  similar  ones  obtained  by  reversing  the  direction  of  the  external  Held  (not 
shown  here)  we  note  there  is  slightly  more  structure  in  the  negative  field  than  in  zero  field 
(not  shown  here)  or  the  positive  field  (Figure  5).  The  structure  is  most  prominent  on  the 
metal  side  because  the  O-H  bond  is  longer  than  PC-0  bond  and  so  the  proton-proton  im¬ 
age  interaction  being  larger  pulls  the  water  molecule  closer  to  the  metal. 

In  the  PC_ST2  density  profile  in  Figure  5  a  distinct  shoulder  occurs  at  about  0.7  nm 
clearly  inside  the  repulsive  region  of  the  9-3  potential.  This  is  permitted  because  the  wall 
potential  acts  the  center  of  the  water  molecule  not  the  component  atoms.  That  the  system 
is  strongly  polarized  as  shown  by  the  z  component  of  the  dipole  density  pp(z).,  and  the 
microscopic  charge  density  Pq(z).  There  is  a  very  pronounced  oscillation  at  the  surface 
in  both  profiles.  Given  the  existence  of  oscillations  in  pwat(z)  these  are  not  unexpected. 
What  is  not  clear  is  whether  they  result  from  the  superposition  of  two  wall  effects.  This 
has  been  checked  with  larger  N  simulations  where  they  occur  to  the  same  degree  implying 
the  oscillation  is  intrinsic  to  one  surface.  The  microscopic  charge  density  has  a  small 
value  in  the  interior  region  of  the  film  and  a  large  oscillation  centered  at  approximately 
0.7  nm.  The  dipole  polarization  shows  large  negative  deviations  at  both  surfaces,  and  a 
uniform  positive  value  across  the  interior.  The  behaviour  of  the  dipole  and  charge  den¬ 
sities  is  that  consistent  with  a  high  dielectric  response  material. 

All  the  simple  rigid  water  models  display  dielectric  saturation  at  high  fields  because  reo¬ 
rientation  is  the  possible  relaxation  process,  and  when  the  molecules  are  all  highly  oriented 
the  dielectric  is  saturated.  In  the  fields  in  the  calculations  used  by  Brodsky  et  al22,  23 
water  was  highly  oriented  and  close  to  saturation,  since  their  system  like  the  one  studied 
here  polarizes  by  molecular  reorientation.  We  have  also  repeated  our  calculations  in  fields 
up  to  1.5  times  stronger  and  find  that  the  interior  polarization  increases  linearly,  implying 
that  dielectric  saturation  has  not  occurred.  However  the  response  of  molecules  near  the 
surface  is  different  and  the  dipolar  density  did  not  increase  linearly  with  the  external  field. 
This  is  consistent  with  the  presence  of  more  highly  oriented  and  more  densely  packed 
water  molecules  next  to  the  surfaces. 

The  question  of  dielectric  response  of  thin  films  and  near  surfaces  in  larger  systems  is  of 
great  interest  If  we  assume  linear  response  of  the  electric  polarization  to  an  external  field 
then  the  dielectric  constant  is  given  by 

T-'-SSrJ 


[14] 


In  the  integrand  Pt  is  the  dipole  density  parallel  to  the  surface  normal  (z)  averaged  over 
the  xy  plane.  The  electrode  charge  is  qT,  and  the  integration  limit  is  bounded  z  <  L/2. 
The  central  portion  of  the  dipole  density  is  linear  in  the  applied  field,  and  appears  to 
provide  a  convenient  'quick'  and  unambiguous  way  to  calculate  the  dielectric  'constant' 
£n.  However  in  practise  this  is  not  easily  realized  because  to  check  convergence  one  must 
evaluate  the  integral  for  very  large  system  to  avoid  interference  from  the  walls. 
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Figure  5.  Water  without  any  ions  between  charged  plates  with  qr=q*c  =  + 1.  Component  density  profile 
plots  for  158  ST2  water  molecules  adjacent  to  charged  metal  electrode  on  right  side.  The  total  electric 
charge  density  on  the  electrode  is  one  positive  electron  on  the  simulation  plate  of  area  L2,  where  L  =  1.862 
nm.  0.29  e/(nm)2  The  vacuum  electric  field  is  approximately  5  GV/m. 

FLUORIDE  SOLUTIONS.  STRONG  ELECTRIC  FIELD  EFFECTS. 

In  this  section  we  discuss  several  simulations  performed  with  one  and  two  fluorides  in  the 
cell  (effective  ionic  concentrations  of  0.35  M  and  0.70  M  respectively)  containing  157  and 
156  water  molecules  respectively.  We  will  show  that  at  the  higher  concentration  there  is 


more  dense  oriented  layer  of  surface  water  and  that  the  inner  surface  electric  field  drives 
the  effect. 

Figure  6  shows  some  representative  results  for  one  fluoride  ion  in  the  form  of  density 
profiles  for  the  fluoride  ion  and  the  atomic  components  of  water.  In  this  calculation 
qr  =  <?<*  =  +  1  and  =  0.  The  surface  field  arises  from  the  negatively  charged  fluoride  ion 
and  its  image.  Because  the  ion  is  distributed  uniformly  in  the  xy  plane  the  surface  electric 
field  is  similar  to  that  inside  a  capacitor  (with  a  1.4  nm  gap).  Features  in  the  water  density 
profiles  near  the  metal  surface  resemble  those  already  seen  for  water  without  ions  in 
Figure  5. 


Figure  6.  Adsorption  of  fluoride  anions.  Density  profiles  for  one  fluoride  F  and  157  ST2  water  molecules 
between  a  metal  electrode  and  the  dielectric  boundary.  Charge  on  the  metal  qr  -  q*  =  +  1 .  Image  plane  at 
z  =  0.931  nm  .  Repulsive  portion  of  the  wall  potentials  begin  at  Izl  =  0.682  nm.  Note  oriented  water  layer 
forming  near  the  metal  (riis)  with  point  charge  (PC  top  curve)  pointing  at  metal. 


Features  on  the  left  at  the  dielectric  boundary  do  not  compare  well,  in  fact  they  resemble 
water  in  zero  field  (not  shown  here).  These  results  show  that  in  the  single  ion  case  there 
is  no  field  at  the  dielectric  and  a  high  field  at  the  metal  surface.  This  indicates  that  the 


inner  field  felt  by  surface  waters  is  similar  to  the  externally  applied  field  in  the  absence 
of  ions.  Closer  comparison  shows  that  at  the  metal  side  the  water  peaks  are  uniformly 
stronger.  In  particular  the  PC_ST2  density  profile  shows  a  distinct  peak  near  0.75  nm,  to 
be  compared  with  the  shoulder  at  0.7  nm  in  Figure  5.  The  extra  height  and  structure  may 
be  due  to  water  in  the  solvation  sheath  of  the  fluoride  ion  at  its  position  of  closest  ap¬ 
proach.  Recall  that  the  fluoride  ion  is  strongly  hydrated  and  shows  no  propensity  to 
contact  adsorb.  Fluoride  remains  hydrated  even  at  the  point  of  closest  approach  (ca.  0.55 
nm)  to  the  metal.  However  at  this  distance  the  solvation  shell  is  splayed  against  the 
electrode.  The  distribution  is  broad  (-0.4  nm  to  0.5  nm,  peak  near  0.2  nm)  covering  almost 
all  accessible  configuration  space  consistent  with  the  ion  being  strongly  hydrated  and  a 
permanent  resident  of  the  diffuse  layer.  Note  also  that  the  fluoride  distribution  is  diffuse 
across  the  entire  water  film,  and  there  being  no  hint  of  a  sharp  plane  or  barrier  (equivalent 
to  an  outer  Helmholtz  plane  OHP  in  the  traditional  model)  across  which  the  ion  is  pre¬ 
vented  from  passing.  In  these  simulations  the  fluoride  behaves  more  like  a  strongly  hy¬ 
drated  ion  in  the  Grahame  model  of  the  electric  double  layer.  Figure  7  shows  density 
profiles  for  the  case  of  two  fluoride  ions  in  the  simulation  cell  with  156  water  molecules. 
The  total  microscopic  charge  density  is  also  shown.  Charge  on  the  metal  qr  =  </«,  =  +  2. 
The  density  profile  for  two  fluorides  resembles  that  of  a  single  fluoride  except  it  is  shifted 
further  from  the  metal  (range  is  from  -0.5  to  0.5  nm  with  the  peak  near  0.0  nm).  The 
electric  Field  due  to  metal  charge  and  ions  is  zero  near  the  dielectric  boundary  as  seen  by 
comparing  water  component  profiles  with  previous  two  figures.  The  inner  surface  field 
has  a  much  higher  value  as  judged  by  sharper  more  intense  peaks  near  the  metal  surface. 
The  sharp  water  peak  near  0.65  nm  has  twice  the  interior  density  value.  Some  of  this 
water  has  one  PC  pointing  at  the  electrode.  Also  note  that  the  fluoride  density  profile  is 
shifted  away  from  the  metal  surface,  consistent  with  the  notion  that  the  compacted  layer 
of  surface  water  hinders  the  approach  of  strongly  hydrated  ions.  This  effect  is  consistent 
with  the  concept  of  die  OHP  being  approximately  two  water  molecules  distant  from  the 
electrode.  The  position  of  the  first  H_ST2  peak,  and  the  near  surface  oscillation  in  mi¬ 
croscopic  charge  confirm  the  overall  picture  as  described. 


Figure  7.  Adsorption  of  anions  continued.  Density  profiles  for  two  fluoride  F  and  156  ST2  water  mole¬ 
cules  between  a  metal  electrode  and  the  dielectric  boundary.  Charge  on  the  metal  qr~<hm  =  +  2.  Image 
plane  at  z  =  0.931  nm.  Repulsive  portion  of  the  wall  potentials  begin  at  Izl  =  0.682  nm.  Note  the  high 
density  oriented  water  peak  near  0.65  nm  and  the  shift  in  the  fluoride  distribution  away  from  the  metal 
compared  to  Figure  6. 


Next  we  show  the  importance  of  the  long  range  interfacial  dipole  field  by  repeating  the 
calculation  for  one  fluoride  shown  in  Figure  6  with  an  uncompensated  field  equivalent  to 
a  charge  of  one  positive  electron  on  the  plate  area  L2.  The  charge  on  the  metal  is  given 
by  qr  =  qim  + qm:  =  +  2,  with  qim  =  q*c  =  +  1.  Figure  8  shows  the  density  profiles  for  one 
fluoride  F  and  157  ST2  water  molecules  between  a  charged  metal  electrode  and  a 
dielectric  boundary.  Note  the  similarities  in  the  density  profiles  for  water  and  its  com¬ 
ponents  in  the  range  z  >  0  nm  in  Figures  7  and  8.  In  this  region  the  density  profile  almost 
superimpose.  Since  we  have  shown  that  the  water  profiles  are  sensitive  to  the  electric  field 
it  suggests  that  die  inner  fields  are  almost  the  same  in  the  two  cases.  Furthermore  these 
calculations  show  that  it  is  the  field  in  the  inner  layer  and  not  steric  effects  due  to  ions 
that  drives  die  structural  changes. 
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Figure  8.  Adsorption  of  anions  continued.  Surface  electric  field  effects.  The  charge  on  the  metal  is  given 
by  <tr  =  <hm  +  q*c  -  +  2,  with  qim  =  quc  =  +  1 .  Density  profiles  for  one  fluoride  F  and  157  ST2  water  molecules 
between  a  charged  metal  electrode  and  a  dielectric  boundary.  Note  similarity  of  water  and  component 
density  profiles  in  Figure  7  for  z  >  0  nm.  Image  plane  at  z  =  0.931  nm.  Repulsive  portion  of  the  wall 
potentials  begin  at  Izl  =  0.682  nm. 


ELECTRIC  FIELD  EFFECTS.  CASE  STUDY  OF  SODIUM  FLUORIDE 

In  this  section  we  consider  the  effect  of  interchanging  anions  and  cations  with  similar 
strengths  of  hydration.  We  have  chosen  NaF  solutions  to  explore  this  case,  and  have 
chosen  to  work  also  with  an  uncompensated  field  (quc  =  ±1)  to  mimic  higher  effective 
concentrations  than  actually  used  in  the  calculation. 

Figure  9  shows  the  density  profiles  for  a  simulation  with  one  fluoride  F,  two  sodium  ions 
Na+  and  155  ST2  water  molecules  between  a  charged  metal  electrode  and  a  dielectric 
boundary.  Hie  charge  on  the  metal  was  qr = <?<#,  +  qIK  =  -2,  with  q*  =  q*  =  - 1.  Note  that 
even  though  the  field  has  opposite  sign  there  is  a  superficial  similarity  of  water  density 
profile  with  those  shown  in  Figure  7  and  8. 


« 


Figure  9.  Approximate  equivalence  interchanged  anions  and  cations,  part  1.  The  charge  on  the  metal  was 
<tr  -  <bm  +  qw  =  -  2,  with  qim  =  qKC  =  -l.  Note  similarity  of  water  and  component  density  profiles  in  Figure 
7.  Density  profiles  for  one  fluoride  F,  two  sodium  ions  Na+  and  155  ST2  water  molecules  between  a 
charged  metal  electrode  and  a  dielectric  boundary.  Image  plane  at  z  =  0.931  nm  .  Repulsive  portion  of  the 
wall  potentials  begin  at  Izl  =  0.682  nm. 


Figure  10  displays  the  profiles  for  two  fluoride  F,  one  sodium  ions  Na+  and  155  ST2 
water  molecules  The  charge  on  the  metal  was  now  positive  <?r  =  <?«.  +  <7kc  =  2,  with 
<lim  -  Quc  =  +  1.  Note  similarity  of  water  and  component  density  profiles  in  Figure  7.  Note 
also  that  the  density  profiles  for  Na  and  F  in  Figures  9  and  10  when  interchanged  are  al¬ 
most  the  same.  These  calculation  extend  to  mixed  electrolytes  the  conclusion  reached  in 
the  last  section  that  die  water  structure  in  the  inner  layer  is  the  result  of  surface  electric 
fields. 
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Figure  10.  Approximate  equivalence  interchanged  anions  and  cations,  pan  2.  Note  similarity  of  water  and 
component  density  profiles  in  Figure  7.  The  charge  on  the  metal  was  positive  qr  -qim  +  q*  =  2,  with 
q*  =  qm:  =  +  1.  Density  profiles  for  two  fluoride  F,  one  sodium  ions  Na+  and  155  ST2  water  molecules 
between  a  charged  metal  electrode  and  a  dielectric  boundary.  Image  plane  at  z  =  0.931  nm  .  Repulsive 
portion  of  the  wall  potentials  begin  at  lz<  =  0.682  nm. 


SMALL  CATION  COADSORBED  WITH  IODIDE  IONS. 

In  this  set  of  simulations  we  explore  another  important  aspect  of  adsorption  on  metal 
electrodes,  namely  the  ability  of  strong  contact  adsorbers  like  iodide  ions  I'  to  adsorb  on 
positively  charged  electrodes  in  sufficient  excess  to  change  the  sign  of  the  charge  at  the 
interface  as  observed  by  an  ion  in  the  diffuse  layer.  In  this  case  cations  are  attracted  out 
of  the  diffuse  layer  region  to  compensate  the  excess  negative  ion  charge  at  the  interface. 
The  lithium  ion  was  chosen  as  cation.  We  have  performed  these  calculations  without  and 
with  uncompensated  charge. 

In  the  first  case  the  charge  on  the  metal  was  qr=qim  =  +  1,  and  in  the  second  case 
qr  =  qlm  +  qm:  =  2.  Figure  1 1  displays  the  density  profiles  for  all  components  of  the  system 
and  the  microscopic  charge  density.  The  simulation  time  was  1000  ps  with  die  first  100 
ps  discarded  for  equilibration,  and  then  configurations  were  stored  every  1.0  ps. 


Note  that  both  iodide  ions  adsorbed  in  one  sharply  peaked  distribution  mostly  inside  the 
onset  of  the  repulsive  wall  region,  but  with  a  small  tail  out  to  smaller  z  positions.  The 
iodide  distribution  in  zero  field  is  characteristically  different  from  the  case  when  an  at¬ 
tractive  field  exists  there  being  a  longer  tail  into  the  electrolyte  indicating  a  more  weakly 
bound  state.  The  single  Li  ion  occupies  a  (possibly  weakly  bimodal)  diffuse-like  distrib¬ 
ution  between  -0.6  nm  and  0.6  nm.  The  water  profiles  hint  that  the  inner  suit  .ce  electric 
field  across  the  first  water  layer  is  weak  because  the  iodide  charge  is  so  close  to  the  metal, 
i.e.,  the  waters  are  outside  the  main  part  of  the  capacitor  is  this  model.  Additionally  the 
large  size  of  the  iodides  tends  to  exclude  some  water  from  the  surface. 
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Figure  11.  Coadsorption  of  ions.  The  charge  on  the  metal  was  qr  =  <ftm-  +  1.  Interface  becomes  effectively 
negatively  charged  after  the  adsorption  of  two  iodide  ions,  thereby  attracting  the  positive  Li  ion.  Density 
profiles  for  two  I\  one  lithium  cation  Li*  and  155  ST2  water  molecules  using  the  immersed  electrode  model. 


In  the  second  case  the  charge  on  the  metal  was  <?r = <7<*  +  =  +  2  due  to  the  presence  of 
the  extra  charge  =  +1.  The  net  charge  on  the  electrode  repells  the  lithium  ion.  Figure 
12  displays  the  density  profiles  for  all  components  of  the  system  and  the  charge  density. 
The  simulation  time  was  1000  ps  with  the  first  100  ps  discarded  for  equilibration,  con¬ 
figurations  were  stored  every  0.S  ps.  Note  that  both  iodide  ions  were  adsorbed  in  a  sharply 


peaked  distribution  that  resembles  the  single  adsorbed  iodide  distribution  calculated  by 
Glosli  and  Philpott16.  The  single  Li  ion  occupies  diffuse-like  distribution  between  -0.6 
nm  and  0.3  nm  and  peaks  around  -0.1  nm.  The  region  between  0.0  and  0.3  nm  defines 
a  region  of  reduced  probability  of  penetrating  into  inner  layer  and  contacting  either  the 
wall  or  the  adsorbed  iodide  ions. 

It  is  noticeable  that  the  water  structure  near  the  metal  surface  though  less  well  defined  than 
in  the  case  of  adsorbed  fluoride  (see  Figures  7,  8  or  10)  shows  more  structure  than  in 
Figure  11.  Overall  the  structure  resembles  water  in  the  field  of  a  positively  charged 
electrode  as  shown  in  Figure  5.  Some  artifacts  are  expected  due  to  some  surface  PC 
reaching  the  high  field  region  between  the  iodide  layer  and  the  metal. 

A  final  word  of  caution.  The  statistics  for  lithium  ions  in  the  simulations  described  in  this 
section  are  not  as  good  as  in  previous  calculations.  The  presence  of  two  iodides  on  the 
surface  creates  a  rough  surface  and  the  calculations  should  be  run  several  nanoseconds  to 
permit  the  positive  ion  to  explore  all  configuration  space. 


Figure  12.  Coadsorption  of  ions  continued.  Density  profiles  for  two  T,  one  Li+  and  1SS  ST2  water  mole¬ 
cules  next  to  immersed  electrode.  The  charge  on  the  metal  was  qr  =  q*  +  q*  =  +  2. 


CONCLUSIONS 


In  this  paper  we  have  shown  how  a  simple  model  suffices  to  mimic  many  phenomena 
familiar  from  experiments  on  electric  double  layers  at  the  electrolyte-metal  interface. 
There  was  emphasis  on  electrostatic  interactions  and  understanding  the  role  played  by  the 
inner  surface  field  in  driving  structural  changes  in  the  surface  water  layer.  This  is  ex¬ 
pected  to  be  important  for  all  polar  systems.  A  key  element  of  the  calculations  was  the 
use  of  the  fast  multipole  method  to  accurately  and  efficiently  calculate  coulomb  inter¬ 
actions  so  that  macroscopic  electric  fields  were  computed  correctly.  Among  the  phe¬ 
nomena  studied  were:  an  oriented  boundary  layer  of  water  at  the  electrode  when  it  is 
charged,  penetration  of  nominally  diffuse  layer  species  like  hydrated  fluoride  into  the  inner 
layer,  and  attraction  of  cations  to  a  positively  charged  electrode  induced  by  contact 
adsorption  of  large  ions  like  iodide. 

Finally  we  point  out  that  the  techniques  described  in  this  report  with  minor  modifications 
can  be  extended  to  more  complex  systems,  for  example:  microelectrodes,  polymer  coated 
electrodes,  biological  membranes,  globular  protein  surfaces,  and  clay  surfaces. 
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